function [y,tau_p,A_p,a_p] = MChannelchirp(y_cha_in,MULTIPATH_NUM,fs,BANDWIDTH)
%% channel delay and attenuation parameters
% c= 1500;                % audio propagation speed in m/s, in water environment
% Nzp=100000;
% Offset =0;
% att_rate = 20/Nzp;       	%  channel attenuates 20 dB within guard time
% path_interval = 0.0015;     % inter-arrival time of paths is set to 1 ms
% path_interval = 1e-3;
% Doppler related parameters
% v_mean = 0;%1.5; 0.1            % 1.50;      % initialize the moving speed here
% v_std = 0;%1;   0.0001           % speed in m/s
% a_mean = v_mean/c;
% a_std = v_std/c;
    
tau_p = [0.0016 0.0042 0.0091  0.0175	0.0276];%EVA model
tau_p=tau_p-tau_p(1);
% generate random channel taps
pdp = [-13.6345	-11.4405	-12.0822	-17.7120	-18.5388];
pow_prof = 10.^(pdp/10);
pow_prof = pow_prof/sum(pow_prof);%normalization of power delay profile
A_p = sqrt(pow_prof).*(sqrt(1/2) * (randn(1,MULTIPATH_NUM)+1i*randn(1,MULTIPATH_NUM)));%channel coef. for each path

% generate Doppler scale factors
%a_p = a_mean*ones(1,MULTIPATH_NUM);
% a_p = a_mean + a_std*randn(1,MULTIPATH_NUM);
% a_p = sort(a_p,'descend');
a_p = [0.003703,	0.003333,	0.002962	,0.003333	,0.003703];
% A_p = ones(1,MULTIPATH_NUM);
% tau_p = zeros(1,MULTIPATH_NUM);
% tau_p  = [0,1e-3];
%% pass through the multipath UWA channel
Npa = MULTIPATH_NUM;
y = zeros(1,1);
% 原始信号长度
y_num = length(y_cha_in);
% 原始信号采样
t0 = 0:1/fs:(y_num-1)/fs;
for p = 1:Npa
    % 受多普勒因子影响后的信号长度
    alfmnum = roundn(y_num/(a_p(p)+1),0);     %用四舍五入函数实现取整操作，这一步是公式提出了a_p+1后信号则缩短了a_p+1倍的长度
    % 新的延时
    tau_new(p) = tau_p(p)/(1+a_p(p));
    % 设计出新的路径信号的时间
    t1=0:1/fs*(y_num/alfmnum):(y_num-1)/fs;%缩短过后的信号的长度
    % 每条径的信号不同的压缩或拓展
    y_multiple = interp1(t0,y_cha_in,t1,'spline');
     % stau = zeros(1,roundn((tau_new(p)-tau_new(1))/(a_p(p)+1)*fs,0));%时延对应到仿真中有多少个点,减去tau_new(1)的目的是为了准确的进行同步
    stau = zeros(1,round(tau_new(p)/(a_p(p)+1)*fs));%如果前面减去了tau_p(1）这里就不用减，
    y_multiple = [stau,y_multiple];
    if length(y_multiple)-length(y)>0
        y =[y,zeros(1,length(y_multiple)-length(y))];
    end
    if length(y_multiple)-length(y)<0
        y_multiple=[y_multiple,zeros(1,length(y)-length(y_multiple))];
    end
    y = y + A_p(p)*y_multiple;
end
end